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I.  SUMMARY 


Four  major  research  tasks  were  performed  during  the  year.  Several  prototype 
liquid  atomization  experimental  designs  were  developed  and  tested,  and  the  pulse 
visualization  system  was  demonstrated.  A  nonlinear  theory  for  liquid  sheet  deforma¬ 
tion  has  been  developed  that  predicts  vortical  rollup  of  the  liquid-gas  interface  and 
allows  for  a  rough  estimate  of  early  ligament  dimensions  in  the  atomization  process. 
The  turbulent  reactive  flow  study  identifies  that  the  mixing  and  entrainment  rates 
are  identical  on  the  high-speed  and  low-speed  sides  of  a  constant-density  mixing 
layer.  The  formation,  merging,  and  pairing  of  vortical  structures  is  also  predicted. 
Variable  property  effects  have  been  included  in  the  vaporizing  droplet  theory  pro¬ 
viding  significant  modifications  on  drag  coefficients.  Nusselt  number,  and  Sherwood 
numbers.  Modified  correlations  for  higher  transfer  numbers  are  being  developed  and 
interacting  droplets  are  considered. 


II.  INTRODUCTION 

The  AFOSR  Program  on  Fundamental  Studies  of  Spray  Combustion  and  Turbu¬ 
lent  Combustion  has  centered  on  four  major  tasks  during  the  reporting  period  from 
1  November  1987  to  31  October  1988.  Significant  progress  has  occurred  on  each 
task.  The  major  objectives  and  the  major  accomplishments  are  summarized  later  in 
this  section.  Details  are  provided  in  the  next  section  entitled  Research  Status.  Of 
course,  even  further  information  is  contained  in  the  publications  listed  in  Section  IV. 

The  three  theoretical  tasks  have  each  yielded  some  very  noteworthy  results  and 
collectively  have  resulted  in  ten  publications  and  four  presentations  at  conferences 
and  seminars  appearing  in  the  last  year  ending  October  31,  1988.  See  Sections  IV 
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and  V.  The  experimental  task  has  resulted  in  a  functioning  apparatus  that  is  yield¬ 
ing  interesting  results.  Further  publications  in  each  of  the  task  areas  are  expected. 
The  program  staff  has  remained  reasonably  stable.  Professors  Samuelsen  an  d  Sirig- 
nano,  Drs.  Rangel  and  Sowa,  and  Messrs.  Chiang,  Miralles- Wilhelm,  Stapper.  and 
Crum  continue  with  the  program.  Dr.  Raju  has  left  for  a  position  elsewhere.  So 
degrees  were  granted  during  the  reporting  period.  The  professional  personnel  are 
listed  in  Section  VI. 

A  summary  of  the  objectives  and  major  findings  and  accomplishments  during 
the  past  year  follows  for  each  task. 

A.  ATOMIZATION  EXPERIMENT 

The  goal  of  the  experimental  program  is  to  observe  and  quantify  the  deforma¬ 
tion  of  a  liquid  sheet  between  the  point  of  injection  and  the  point  of  breakup.  Of 
specific  interest  are  the  wave  length  and  wave  amplitude  of  the  deformed  sheet.  The 
objective  of  the  current  year  has  been  to  establish  and  demonstrate  an  experimen¬ 
tal  fixture  that  produces  the  desired  condition,  and  to  establish  and  demonstrate  a 
visualization  system  for  conducting  the  measurements. 

During  the  year,  a  prototype  test  fixture  (labeled  Prototype  II)  was  fabricated 
based  on  design  changes  implemented  to  reduce  the  pressure  drop  in  the  nozzle  in¬ 
lets  and  to  improve  the  machinability  of  the  fixture.  The  prototype  was  tested  and 
found  to  produce  a  sheet  of  higher  quality  than  that  obtained  using  the  Prototype  I 
fixture.  The  quality  of  the  sheet,  however,  was  not  sufficient  to  meet  the  goals  of  the 
program.  Longitudinal  striations  in  the  formed  sheet  continued  to  persist  although 
a  substantial  improvement  was  noted  over  Prototype  I.  Three  scenarios  have  been 
identified  to  obtain  the  quality  of  sheet  desired:  (1)  Reduce  surface  roughness  within 
the  nozzle,  (2)  increase  the  liquid  flow  velocity,  and  (3)  improve  the  internal  nozzle 


hydrodynamics. 


By  the  end  of  the  current  reporting  period,  the  first  scenario  was  addressed  with 
favorable  results.  Teflon  was  sprayed  on  the  internal  surfaces  of  Prototype  II.  and 
a  third  experimental  fixture  (labeled  Prototype  III)  was  machined  from  INCONEL 
to  improve  material  strength  and  resistance  to  corrosion,  and  the  surfaces  were  ma¬ 
chined  to  minimize  surface  roughness.  Improvements  in  the  quality  of  the  Liquid 
sheet  were  clearly  evident.  At  the  end  of  the  reporting  period,  steps  were  initiated 
to  address  the  second  and  third  scenarios,  both  of  which  are  likely  contributors  to 
the  sheet  quality. 

The  pulse  visualization  system  was  demonstrated  using  a  non-  TV-optimized 
imaging  lens.  The  system  will  become  fully  operational  when  a  custom  made  TV 
lens  package  arrives  early  in  the  next  reporting  period.  With  the  TV  lens  package 
in  place,  the  system  will  be  used  both  in  a  continuous-filming  and  a  single-shot 
mode,  the  latter  of  which  will  be  used  to  perform  the  computer  image  enhancement 
routines  and  sizing  calculations  for  the  sprayed  sheet. 

B.  ATOMIZATION  1.HEQRY 

It  is  well  known  that  the  efficiency  of  a  liquid-fueled  combustor  depends  criti 
callv  on  the  efficiency  of  the  vaporization  process  which  is  strongly  dependent  on  the 
droplet  size  distribution  resulting  from  the  atomization  mechanism.  Deterministic 
models  of  atomization  processes  are  still  in  their  infancy,  even  for  the  simplest  con¬ 
figurations.  Some  empirical  models  based  on  a  few  more  or  less  ad  hoc  assumptions 
exist,  but  the  simulation  of  the  evolution  of  a  large  liquid  mass  going  through  a 
nonlinear  deformation  process  and  breakup  eventually  resulting  in  a  spray  is  still 
to  come.  Our  efforts  in  this  area  are  devoted  to  the  study  of  the  nonlinear  shear 
distortion  of  a  planar  liquid  sheet  interacting  with  a  stream  of  air  on  each  side. 
The  effect  of  surface  tension  as  well  as  that  of  a  density  dicontinuity  are  included. 
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Linear  theory  results  have  been  generalized  to  include  the  finite-sheet  case.  The 
critical  wavenumber  is  shown  to  be  weakly  dependent  on  the  sheet  thickness  when 
the  density  ratio  is  very  small,  as  in  air-liquid  systems.  The  numerical  calculations 
show  qualitative  differences  as  the  sheet  thickness  decreases  or  as  surface  tension 
increases,  mainly  a  disappearance  of  the  rollup  features  encountered  in  the  single¬ 
sheet  case.  Some  dimensions  of  the  early  ligaments  in  the  breakup  process  can  be 
estimated. 


C.  TURBULENT  REACTIVE  FLOWS 

The  two-way  interaction  between  the  fluid  dynamics  and  the  rate  of  energy  re¬ 
lease  plays  a  fundamental  role  in  the  evolution  of  turbulent  reactive  flows.  The 
objective  of  this  study  is  to  gain  a  better  understanding  of  the  transport  processes 
that  occur  in  turbulent  flows,  particularly,  the  interaction  and  relative  importance 
of  the  different  rate-controlling  factors.  Mixing  rate  and  reaction  rate  predictions 
for  turbulent  mixing  layers  are  a  further  objective.  Developments  during  the  re¬ 
ported  year  have  occurred  along  two  directly  related  fronts.  On  one  hand,  the 
species  concentration  and  probability  density  function  distributions  in  an  array  of 
non-interacting  viscous  vortices  have  been  obtained  by  an  approximate  analytical 
approach.  When  comparing  the  results  with  experiments  reported  in  the  literature, 
it  is  found  that  the  lack  of  symmetry  shown  in  the  experiments  is  associated  with  a 
bias  of  the  measurements,  and  not  with  real  differences  in  mixing  rate.  Results  for 
the  mixedness  parameter  indicate  that  the  mixing  process  is  enhanced  by  increasing 
the  vortical  strength,  and  weakly  influenced  by  molecular  diffusi vity. 

The  second  approach  being  developed  relaxes  the  assumption  of  noninteractive 
vortical  structures  by  considering  the  evolution  of  the  flow  field  from  its  initial  con¬ 
figuration  consisting  of  a  surface  of  discontinuity  to  the  formation  and  interaction 
of  large  scale  vortical  structures.  During  the  phase  reported  here,  the  efforts  have 
been  concentrated  on  understanding  the  fluid  mechanics  that  control  the  pairing 
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and  merging  process  and  the  ways  in  which  it  can  be  controlled  by  means  of  appro¬ 
priate  external  excitation.  It  is  shown  that  the  interaction,  grouping,  and  eventual 
merger  of  large-scale  vortical  structures  is  the  result  of  basically  two  distinct  phe¬ 
nomena:  the  shearing  effect  of  the  initial  layer  configuration  which  tends  to  bring 
the  structures  together  so  that  pairing  can  occur,  and  the  spreading  effect  of  viscous 
diffusion  which  is  responsible  for  the  merger  of  the  structures. 

D.  VAPORIZING  DROPLET  CALCULATIONS 

The  general  objective  of  this  study  is  the  detailed  analysis  of  a  vaporizing  droplet 
or  an  interesting  group  of  vaporizing  droplets  in  a  convective  environment.  In  or 
der  to  determine  the  dimensions  and  predict  performance,  stability,  and  pollutant 
emission  of  a  combustor  design,  the  trajectory,  heating  and  vaporizing  history  of  a 
droplet  must  be  accurately  computed.  Hence  this  research  has  been  focusing  on  the 
accurate  investigation  of  the  local  transport  processes  as  well  as  the  overall  param¬ 
eters  describing  the  system  such  as  net  drag  force,  heat  and  mass  transfer  around 
the  droplet.  It  is  also  desirable  to  obtain  the  relationship  between  the  transfer  co¬ 
efficients  and  some  important  fluid  parameters  such  as  heat,  mass  transfer  numbers 
and  Reynolds  number. 

The  unsteady  Navier-Stokes  equations  are  numerically  solved  with  minimum  as¬ 
sumptions.  The  numerical  calculation  is  performed  in  generalized  coordinates  over 
a  computational  domain,  which  is  accordingly  adjusting  to  the  moving  interface 
boundary.  The  effect  of  variable  properties,  which  is  believed  to  be  very  important 
in  a  typical  high  pressure  and  high  temperature  flow  field,  is  the  main  consideration 
during  this  research.  The  results  indicate  that  the  constant  properties  case  over¬ 
predicts  the  drag  coefficient  at  least  by  20%.  The  drag  coefficient  does  not  simply 
increase  with  time  as  a  result  of  the  reduction  of  Reynolds  number.  It  is  noteworthy 
that  the  transfer  number  plays  a  very  important  role  in  the  droplet  transport  phe¬ 
nomena.  For  higher  transfer  numbers,  the  vaporization  rate  is  larger  and  the  drag 
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coefficient  can  be  significantly  reduced.  For  lower  transfer  numbers,  the  boundary 
layer  blowing  effect  is  weaker  and  the  droplet  acts  more  like  as  a  solid  sphere  while 
the  drag  coefficient  is  dominated  by  the  Reynolds  number  reduction.  The  correla¬ 
tion  of  Renksizbulut  and  Haywood  does  not  agree  with  our  results  for  the  range  of 
transfer  numbers  considered  in  the  present  calculation. 

We  have  examined  the  interacting  effects  between  two  vaporizing  droplets  which 
are  moving  in  tandem  with  respect  to  the  free  stream,  in  an  intermediate- Reynolds- 
numbers-flow  over  a  limited  range  of  different  initial  Reynolds  number,  droplet  spac¬ 
ing.  and  droplet  radii  ratio.  A  very  sophisticated  grid  generation  technique  is  em¬ 
ployed  to  accommodate  the  changing  boundary  shapes.  The  basic  solution  procedure 
is  the  same  with  the  single  droplet  case.  The  constant  properties  assumption  (ex¬ 
cept  gas-phase  density)  is  used  to  simplify  the  algorithm  involved  solving  interface 
boundary  conditions  simultaneously.  The  results  show  that  the  droplet  interactions 
are  evident  for  initial  Reynolds  numbers  of  50  to  200  and  for  initial  droplet  spacing  of 
2  to  15  droplet  diameters.  The  droplets  can  either  collide  or  move  apart  depending 
upon  the  initial  conditions.  There  exists  a  critial  ratio  of  initial  droplet  diameters 
a2.o/a\,o  below  which  droplet  collision  becomes  unlikely.  The  momentum,  heat  and 
mass  transfer  processes  are  quite  different  for  two  droplets  depending  upon  the  spac¬ 
ing.  For  a  sufficient  spacing,  both  droplets  behave  as  two  isolated  droplets.  When 
the  droplet  spacing  becomes  too  small,  the  downstream  droplet  is  fully  covered  by 
the  wake  of  the  lead  droplet  and  the  heat  transfer  is  significantly  reduced. 


III.  RESEARCH  STATUS 

A.  ATOMIZATION  EXPERIMENT 


Experimental  Fixture.  During  the  past  year,  the  second  prototype  fixture  (la¬ 
beled  Prototype  II)  was  fabricated,  tested  and  evaluated,  modified,  and  reevaluated. 


The  design  of  this  new  two-dimensional,  water-sheet  nozzle  was  developed  from  in¬ 
formation  obtained  through  the  testing  and  evaluation  of  the  Prototype  I  fixture. 
The  goals  of  the  Prototype  II  design  were  to  reduce  the  pressure  drop  in  the  nozzle 
inlets  and  to  improve  the  machinability  of  the  nozzle.  The  Prototype  II  fixture  was 
machined  from  6061  aluminum  to  facilitate  machining  and  a  higher  quality  surface 
finish  than  the  2024  aluminum  used  in  the  first  prototype.  The  Prototype  II  was 
fabricated  wfith  increased  clearances  within  the  nozzle  to  reduce  the  pressure  drop 
prior  to  the  nozzle  tip.  and  with  added  supports  to  maintain  constant  gap  thickness 
across  the  width  of  the  nozzle. 

The  Prototype  II  fixture  was  tested  using  the  same  flow  settings  that  were  used 
with  the  Prototype  I.  The  initial  results  were  promising,  and  much  improved  over 
the  performance  of  the  Prototype  I.  However,  the  liquid  sheet  was  not  completely 
satisfactory.  Further  testing  of  thePrototype  II  resulted  in  a  gradual  deterioration 
of  the  quality  of  the  sheet.  The  sheet  was  found  to  be  perturbed  by  rough  surface 
effects  within  the  nozzle.  Inspection  of  the  nozzle  revealed  that  the  6061  aluminum 
oxidized  readily  in  the  presence  of  the  water.  The  oxidation  resulted  in  a  substan¬ 
tially  increased  roughness  of  the  surfaces  within  the  nozzle. 

The  Prototype  II  fixture  was  modified  by  (1)  increasing  the  thickness  of  the 
liquid  sheet  to  minimize  sheet  perturbations,  (2)  coating  the  fixture  surfaces  with 
Teflon  spray  to  further  smooth  the  surface  and  prevent  oxidation  of  the  aluminum, 
and  (3)  polishing  the  interior  filming  lip  of  the  fixture  with  sub- micron  lapping  film 
to  reduce  surface  roughness  effects.  The  increase  in  the  thickness  of  the  liquid  sheet 
reduced  the  size  of  the  surface  imperfections  relative  to  the  sheet  thickness.  The 
coating  of  the  surfaces  with  Teflon  created  a  surface  that  not  only  resisted  oxida¬ 
tion,  but  also  had  a  high  contact  angle  with  the  water.  Polishing  the  inner  surfaces 
with  sub-micron  lapping  film  produced  a  mirror-like  finish.  These  modifications  im¬ 
proved  the  quality  of  the  sheet.  However,  it  was  evident  that  the  quality  of  the  sheet 
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was  still  not  sufficient  to  meet  the  goals  of  the  program.  While  Prototype  II  was 
being  modified  and  tested,  a  new  experimental  fixture.  Prototype  III.  was  designed 
and  fabricated  from  INCONEL  in  order  to  eliminate  oxidation  of  the  inner  surfaces. 
Thicker  shims  were  specified  to  increase  the  sheet  thickness  by  a  factor  of  five.  The 
width  of  the  sheet  was  also  increased  by  509c  from  Prototype  II  to  Prototype  III  to 
allow  for  greater  ease  in  taking  measurements. 

The  Prototype  III  was  demonstrated  at  liquid  velocities  ranging  from  3.3  to  7.4 
m/sec,  and  shear  air  velocities  from  0  to  54  m/sec.  The  tests  were  conducted  with 
shear  air  on  both  sides  of  the  sheet,  and  on  only  one  side  of  the  sheet.  The  sheet 
thickness  for  the  tests  was  508  pm  and  the  sheet  width  was  3S.1  mm.  The  charac¬ 
teristic  measurement  obtained  during  the  testing  was  the  distance  from  the  tip  of 
the  nozzle  to  the  point  at  which  droplets  began  breaking  off  from  the  sheet. 

The  results  of  these  measurements,  plotted  versus  the  shear  air  velocity,  are 
given  in  Figures  A.l  and  A. 2.  In  Figure  A.l,  the  distance  to  breakup  is  plotted 
versus  the  air  velocity  at  a  series  of  liquid  flow  rates.  All  cases  appear  to  approach 
an  asymptotic  value  as  the  air  velocity  decreases.  This  is  expected  as  the  breakup 
of  the  sheet  is  dependent  upon  the  relative  velocity  of  the  liquid  sheet  and  the  shear 
air.  As  the  air  velocity  decreases,  the  air-  liquid  system  approaches  zero  relative 
velocity,  at  which  point  there  is  no  shearing  to  set  up  surface  waves  on  the  sheet, 
and  therefore,  no  breakup  of  the  sheet.  In  theory,  as  the  liquid  velocity  is  increased 
above  the  air  velocity,  the  sheet  at  some  point  will  breakup.  The  limit  for  a  given 
liquid  velocity  is  when  the  air  velocity  approaches  zero.  This  was  not  observed 
during  the  Prototype  III  shake-down  tests  possibly  because  the  water  momentum 
was  not  sufficiently  strong  to  overcome  the  nozzle-geometry-  induced  sheet  imper¬ 
fections.  To  evaluate  this  hypothesis,  the  flow  capacity  of  the  system  is  presently 
being  increased  to  allow  for  testing  with  liquid  velocities  up  to  50  m/sec. 
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A  second  observation  from  these  tests  is  that  the  distance  from  the  nozzle  tip 
for  sheet  breakup  increases  as  the  liquid  velocity  increases  for  a  given  air  flow  rate. 
In  Figure  A. 2,  the  distance  from  the  nozzle  tip  to  the  point  that  droplets  break  off 
the  sheet  is  plotted  versus  the  shear  air  velocity  for  a  liquid  velocity  of  7.4  m/sec. 
The  tv  v  cases  are  shown:  shear  air  flow  on  one  side  of  the  liquid  sheet,  and  shear 
air  flow  on  both  sides  of  the  sheet.  In  both  instances,  the  distance  from  the  nozzle 
tip  to  the  breaking  off  of  the  droplets  decreases  with  increasing  air  velocity.  An 
interesting  feature  of  Figure  A. 2  is  that  the  distance  to  breakup  is  constant  for  the 
same  total  air  flow.  In  other  words,  in  the  two-sided  flow  case,  the  air  velocity  is 
only  half  the  air  velocity  required  in  the  single-sided  case  for  the  same  distance  to 
sheet  breakup.  However,  since  the  breakup  described  in  the  model  is  based  upon 
the  relative  velocity  of  the  liquid  and  air.  this  breakup  at  low  liquid  velocity  may  be 
due  to  another  cause.  That  is.  the  momentum  of  the  air  component  normal  to  the 
liquid  sheet  may  be  the  dominating  factor  influencing  the  breakup,  in  which  case 
it  would  be  expected  that  the  distance  to  breakup  would  be  the  same  for  the  same 
total  air  flow  rate. 

Future  Direction.  In  our  efforts  to  improve  the  quality  of  the  liquid  sheet,  we 
have  come  across  some  new  ideas  from  fixtures  used  in  large  scale  lasers  that  flow  a 
two-dimensional  sheet  of  dye  solution  through  the  laser  cavity.  In  laser  applications, 
the  flow  undergoes  constant  acceleration  until  the  nozzle  exit.  This  is  illustrated  in 
Figure  A. 3.  It  is  anticipated  that  this  modification  will  provide  liquid  sheet  qualities 
better  than  those  displayed  by  the  earlier  prototypes. 

Pulse  Visualization  System.  The  pulse  visualization  system  has  been  success¬ 
fully  demonstrated  using  a  non-UV-  optimized  imaging  lens  (Canon  18-108mm  1:2.5 
Zoom  lens  with  a  6X  zoom  ratio).  This  lens,  while  satisfactory  as  a  test  for  the  basic 
system,  attenuates  the  lTV  signal  from  the  nitrogen  laser  and  can  only  be  used  in 
the  visible  light  range.  In  addition,  the  mechanical  shutters  are  not  fast  enough  to 
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stop  the  action  of  the  sprayed  sheet. 


The  Particle  Measuring  System  custom-made  UV  lens  package  will  be  the  per¬ 
manent  imaging  system  used  with  the  pulse  visualization  system  and  is  on  order. 
The  short  pulse  duration  of  the  nitrogen  laser,  coupled  with  the  UV  lens  package, 
is  ideal  for  taking  sharp,  clear  pictures  of  the  spray.  The  laser  will  be  used  pri¬ 
marily  in  the  single-shot  mode  for  most  of  the  data  collection,  so  that  individual 
frames  can  be  manipulated  by  the  computer  through  image  enhancement  routines. 
However,  for  imaging  the  spray  in  real-time,  without  computer  enhancement,  the 
system  can  be  run  up  to  20  Hz,  or  20  frames  per  second.  The  camera  can  also  be 
used  to  film  on  a  continuous  basis.  When  the  camera  is  operating  in  the  continuous 
filming  mode:  however,  it  does  not  update  the  image  with  a  high  enough  frequency 
to  obtain  high-resolution  pictures  of  the  moving  spray. 

B.  ATOMIZATION  THEORY 

The  effort  during  this  period  was  aimed  at  extending  the  mathematical  and  nu¬ 
merical  model  described  in  the  last  annual  report  (see  also  Rangel  and  Sirignano. 
1988)  in  order  to  simulate  the  behavior  of  a  finite-thickness  liquid  sheet  in  contact 
with  a  stream  of  air  on  each  side  of  it.  In  this  context,  we  should  note  that  our 
previous  analysis  of  a  single  interface  represents  the  limit  of  an  infinitely  thick  sheet 
or.  more  realistically,  a  case  in  which  the  sheet  thickness  is  much  larger  than  the 
disturbance  wavelength  and  the  disturbance  amplitude.  Linear  theory  may  be  used 
also  in  the  finite-sheet  case  to  predict  the  initial  growth  of  disturbances  on  the  liquid 
sheet.  Figure  B.l  is  a  schematic  of  two  limiting  situations  that  can  be  investigated. 
We  first  consider  cases  in  which  the  air  streams  on  either  side  of  the  liquid  sheet 
move  with  the  same  velocity.  In  one  case  shown  in  Fig.  B.l,  symmetric  waves  are 
studied.  In  the  other  case,  asymmetric  or,  more  properly,  shifted-symmetric  waves 
are  considered.  Much  of  the  effort  during  the  reported  period  has  concentrated  on 
generalizing  linear  theory  results  to  be  used  as  the  starting  point  for  the  numerical 
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calculations. 


The  situation  examined  at  first  consists  of  a  planar  Liquid  sheet  of  finite  thick¬ 
ness  in  contact,  on  each  side,  with  a  stream  of  air  of  thickness  much  larger  than 
that  of  the  liquid  sheet.  Briefly,  inviscid  linear  theory  reduces  the  problem  to  the 
solution  of  Laplace's  equation  in  each  stream  subjected  to  the  interface  conditions 
of  continuity  of  normal  velocity  and  momentum  balance.  Modal  analysis  is  used  to 
represent  the  evolution  of  either  interface  as: 

r;  =  e'ul  sin(2~j  )  (1) 

Linear  theory  shows  that  the  shifted-symmetric  waves  are  more  unstable  and 
for  such  a  case  it  predicts  a  cutoff  wavenumber  that  divides  the  spectrum  into  a 
region  of  stable  waves  and  a  region  of  exponentially  growing  waves.  This  cutofT 
dimensionless  wavenumber  is  found  to  be: 


1  +  pCOth(27T/i/A) 

where  p  is  the  density  ratio  (air/liquid),  h  is  thickness  of  the  undisturbed  sheet  and 
A  is  the  disturbance  wavelength.  Two  things  are  immediately  verified  from  Eq.  1. 
First,  the  limit  h  — >  oc  yields  U’c  =  1  as  corresponds  to  the  single-interface  case. 
Secondly,  the  limit  h  —*  0  results  in  M,-  — >  0  which  indicates  that  the  stability  of  the 
sheet  increases  as  its  thickness  is  decreased.  It  should  be  realized,  however,  that  this 
conclusion  is  valid  only  in  the  case  in  which  the  air  has  the  same  velocity  on  both 
sides.  If  this  were  not  the  case,  an  infinitely  thin  liquid  sheet  would  still  be  unstable 
since  in  such  a  case  the  two  air  streams  would  approach  the  limit  of  a  vortex  sheet 
as  the  liquid  thickness  approaches  zero.  As  in  the  one-sided  case  (infinitely  thick 
sheet),  linear  theory  shows  that  there  is  a  value  of  the  dimensionless  wavenumber 
that  yields  the  maximum  growth  rate.  This  value  is: 
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(3) 


...  _  2(1  +  p)tanh(27rfr/A) 

°vt  3(tanh(2?r/i/A)  +  p) 

The  limit  of  an  infinitely-thick  liquid  sheet  is  rapidly  approached  as  h  increases  and 

yields  the  density  independent  value  of  2/3. 

Our  non-linear  investigation  is  based,  as  it  was  in  the  single-interface  case,  on 
the  vortex  discretization  method.  It  is  described  in  detail  by  Rangel  and  Sirignano 
( 19S8 )  and  in  the  previous  AFOSR  annual  report.  The  existence  of  a  second  inter¬ 
face  means  that  a  second  row  of  vortices  must  be  considered  located  a  distance  h 
below  the  first  one.  This  two-fold  increase  in  the  number  of  vortex  elements  trans¬ 
lates  into  approximately  a  four-fold  increase  in  computational  time.  Surface-tension 
and  density-ratio  effects  are  included  in  the  form  of  an  evolution  equation  for  the 
vortex  strength  of  each  interface.  Other  details  of  the  procedure  are  basically  the 
same  ones  described  in  the  previous  report  and  need  not  be  repeated  here. 

Some  results  are  shown  in  Figs.  B.2.  through  B.6.  These  figures  show  the  evo¬ 
lution  of  disturbances  for  different  configurations  of  interest.  Figure  B.2  shows  a 
situation  that  approaches  the  infinitely-thick-sheet  behavior.  Here  the  disturbance 
wavelength  is  one  half  of  the  liquid  sheet  thickness  (note  that  the  x  and  y  scales  are 
different)  and  each  interface  evolves  almost  independently  of  the  other.  Therefore 
the  rollup  proceeds  as  in  that  case.  When  the  liquid  sheet  thickness  is  reduced  so 
that  the  disturbance  wavelength  is  twice  the  liquid  thickness  (Fig.  B.3).  the  effect 
of  the  other  interface  begins  to  become  noticeable.  In  the  limit  of  zero  thickness 
the  two  air  streams  axe  stable  because  they  have  the  same  velocity.  When  surface 
tension  is  present,  the  rollup  features  are  gradually  eliminated.  At  the  most  unsta¬ 
ble  case,  the  liquid  sheet  evolves  as  illustrated  in  Fig.  B.4,  where  it  may  be  argued 
that  two  ligaments  will  be  produced  for  each  wavelength.  The  symmetric  case  is 
illustrated  in  Figs.  B.5  and  B.6  for  two  different  sheet  thicknesses.  The  rollup  fea¬ 
tures  again  disappear  as  the  thickness  is  decreased  and  only  one  ligament  is  formed 
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for  each  wavelength  of  the  disturbance.  Thus,  the  symmetric  mode  appears  to  yield 
larger  liquid  ligaments  and  presumably  larger  droplets  than  the  shifted-symmetric 
mode. 

Further  developments  are  expected  to  include  a  parametrization  of  the  surface- 
tension  and  density-ratio  effects  for  the  finite-thickness  liquid  sheet  perturbed  in 
the  symmetric  and  shifted-symmetric  modes.  Extension  to  the  case  of  the  two  air 
streams  moving  with  different  velocities  should  follow.  By  studying  the  non-linear 
evolution  of  these  disturbances,  it  is  hoped  that  some  inferences  will  be  made  as  to 
what  the  size  of  the  characteristic  ligaments  is  and  what  effects  the  various  different 
parameters  have  on  this  size. 


C.  TURBULENT  REACTIVE  FLOWS 

The  analysis  of  concentration  distributions  in  a  non-reacting  shear  layer  contin¬ 
ued  during  this  reporting  period.  In  this  new  effort,  the  investigation  focused  on  the 
problem  of  the  interaction  of  a  vortical  flow  with  molecular  species  diffusion,  consid¬ 
ering  the  situation  where  an  infinite  row  of  two  dimensional  vortices  is  superposed 
onto  an  interface  initially  separating  two  species.  Under  the  action  of  the  vortical 
field,  layers  of  both  fluids  tend  to  roll  up  as  shown  in  Fig.  C.l. 

The  level  of  resolution  required  to  quantify  the  molecular  diffusion  process  results 
in  highly  time  consuming  computations.  Because  of  this,  an  approximate  analyti¬ 
cal  approach  based  on  the  methodology  of  Marble  (1985),  Karagozian  (1982),  and 
Karagozian  and  Marble  (1986)  is  employed  to  solve  the  temporal  growth  of  the 
mixing  layer  locally.  This  approach  was  recently  used  by  Cetegen  and  Sirignano 
(1987,1988)  to  determine  the  concentration  field  and  probability  density  functions 
using  Taylor’s  frozen  flow  assumption.  This  analysis  has  been  extended  to  con¬ 
struct  the  pdf  in  a  laboratory  frame  of  reference,  providing  ways  for  a  more  realistic 
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comparison  with  results  of  mixing  layer  experiments  reported  in  the  literature  (Di- 
motakis,  1986;  Dimotakis  and  Brown,  1976;  and  Masutani  and  Bowman,  1986). 

The  velocity  field  for  an  infinite  row  of  equidistant  two  dimensional  vortices 
(Lamb,  1945)  was  modified  to  account  for  a  viscous  core  in  cases  of  one,  three,  or 
five  vortices.  Consideration  of  more  than  five  viscous  vortices  proved  unnecessary,  as 
it  is  reported  by  Miralles- Wilhelm  et  al  (1989).  A  nondimensional  vortex  strength, 
represented  by  the  Reynolds  number,  Re  =  Y /2ni/  .  evolves  form  this  analysis  as  an 
important  parameter  governing  the  process. 

An  analytical  solution  of  the  conservation  of  species  equation  for  a  material  ele¬ 
ment  located  instantaneously  at  the  interface  separating  the  two  non-reacting  species 
can  be  written  in  terms  of  an  error  function  of  appropriately  stretched  distance  and 
time  coordinates,  accounting  for  the  strain  history  of  each  material  element.  The 
nondimensionalization  of  this  equation  introduces  the  Schmidt  number  Sc  -  u/D 
as  the  parameter  controlling  the  molecular  diffusion  process.  The  strain  rate  is  com¬ 
puted  locally  from  the  velocity  field  and  the  shape  of  the  interface. 

Figure  C.2  shows  the  distortion  of  the  initially  horizontal  material  line  separat¬ 
ing  the  two  species  at  certain  instants  after  the  vortices  are  initiated.  The  aging 
process  of  the  vortices  is  taken  into  account,  since  the  vortices  are  born  periodically 
in  time.  Notice  that  as  time  passes,  the  roll  up  of  layers  of  the  fluid  streams  is 
augmented. 

Figure  C.3  shows  the  time  evolution  of  the  concentration  profile  for  the  single 
vortex  case  (Re=50;Sc=l),  as  it  would  be  measured  by  a  probe  located  at  different 
heights  above  the  mixing  plane  originally  dividing  the  two  species.  Notice  that  as 
the  height  y  of  the  probe  decreases  approaching  the  mixing  plane,  a  well  mixed 
region  is  encountered  (i.e.  the  vortex  core).  It  is  in  this  region  where  the  diffusion 
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process  occurs  more  efficiently,  as  expected.  On  the  other  hand,  when  the  height 
of  the  probe  increases,  an  increasing  tendency  towards  non-mixing  can  be  observed, 
until  no  mixing  at  all  is  reached  far  away  from  the  vortex  center.  In  a  frame  of 
reference  attached  to  the  vortical  center,  it  has  been  shown  that  the  concentration 
profiles  exhibit  an  antisymmetry  (Cetegen  and  Sirignano,  1987,  1988).  That  an¬ 
tisymmetry  is  not  found  in  Figure  C.3,  for  a  frame  of  reference  that  is  fixed,  as 
in  an  experimental  situation.  It  is  from  these  curves  that  the  probability  density 
function  at  different  y  locations  is  constructed.  Figure  C.4a  shows  the  pdf  for  the 
case  described  above.  As  it  was  noted  before,  the  higher  peaks  correspond  to  the 
unmixed  concentrations  far  away  from  the  vortex  center,  and  finite  probabilities  of 
intermediate  concentrations  exist  in  the  vortex  core  region,  as  it  is  shown  by  the 
other  relative  smaller  peaks. 

The  numerical  results  for  the  frozen  flo%v  approximation  exhibit  a  characteris¬ 
tic  antisymmetry  whereby  pdf(y.C)  =  pdf(-y.  1  -  C),  while  experimental  evidence 
shows  that  the  peaks  that  occur  in  the  high  speed  stream  side  of  the  flow  are  greater 
than  those  on  the  low  speed  stream  side.  For  this  reason,  the  relative  difference  R  be¬ 
tween  the  pdf  at  a  given  (y.C)  pair  and  its  antisymmetrical  counterpart  (-y.  1  -  C) 
is  plotted  in  Fig.C.4b.  The  numerical  calculations  and  the  shape  of  the  graph  reveal 
three  important  facts.  First  of  all,  that  R  can  be  nonzero  and  therefore  antisym¬ 
metry  is  not  found.  Secondly,  the  pdf  peaks  for  the  high  speed  stream  are  not 
always  greater  than  those  for  the  low  speed  stream  (R  can  be  a  negative  quantity). 
The  other  important  issue  showed  by  this  graph  is  that  the  absolute  value  of  R 
is  much  greater  when  R  is  positive.  This  means  that  when  the  peaks  in  the  high 
speed  stream  are  greater  than  those  on  the  low  speed  stream,  the  relative  difference 
between  them  is  much  greater  (of  the  order  of  5  times)  than  that  in  the  opposite 
case  (when  the  low  speed  stream  peaks  are  higher  than  the  high  speed  stream  ones). 
This  last  result  tends  to  agree  with  the  experimental  results  presented  by  Masutani 
and  Bowman  (1986),  while  the  first  results  support  the  idea  that  these  experimental 
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results  (R  always  positive)  are  a  consequence  of  the  measurement  bias,  rather  than 
a  result  of  a  real  difference  in  mixing  rate.  The  high  speed  side  of  the  vortical  struc¬ 
ture  moves  faster  on  average  towards  the  probe  than  the  low  speed  side,  so  that  the 
probe  sees  a  different  concentration  profile  and  a  different  pdf  from  that  seen  in  a 
frame  of  reference  attached  to  the  vortical  structure.  A  theoretical  explanation  of 
the  influence  of  a  Galilean  transformation  on  the  pdf  is  given  by  MiraUes-Wilhelm 
et  al  ( 19S9). 

The  influence  of  Re  and  Sc  on  the  mixing  process  can  be  presented  in  terms  of 
a  local  or  averaged  mixedness  parameter,  given  by  Miralles-Wilhelm  et  al  (1989) 
and  Cetegen  and  Sirignano  (1987),  respectively.  In  both  cases,  an  increase  in  Re. 
enhances  the  mixing  process  by  increasing  the  rolling  up  of  the  fluid  layers  into 
each  other.  Also,  in  both  cases,  a  weak  dependence  on  Sc  is  observed.  This  results, 
presented  by  Miralles-Wilhelm  et  al  (1989)  are  in  agreement  with  those  presented 
by  Cetegen  and  Sirignano  (1987).  In  a  frame  of  reference  attr  cited  to  the  vortical 
structure,  an  antisymmetric  pdf  results,  and  thus  the  mixedness  results  for  +y  and 
-y  collapse  into  one  curve,  as  observed  in  Fig.  C.5.a.  The  presence  of  the  surround¬ 
ing  vortices  tend  to  decrease  the  rate  of  mixing  away  from  the  mixing  plane,  which 
subsequently  diminishes  the  overall  mixing  process.  This  is  shown  in  Fig  C.5.b. 

In  conclusion,  both  an  increase  in  the  Reynolds  number  and  a  decrease  in  the 
Schmidt  number  tend  to  increase  the  probabilities  in  the  mixed  core  region.  That 
is,  the  mixedness  in  the  layer  is  greatly  enhanced  by  increasing  the  vortical  strength, 
and  its  relatively  poorly  influenced  by  molecular  diffusion.  The  pdf  computations 
relaxing  the  frozen  flow  assumption  suggest  that  the  bias  present  in  experimental 
results  is  not  caused  by  differences  in  mixing  rate.  It  is  not  a  mixing  bias,  it  is  a 
measurement  bias.  This  fact  can  be  explained  theoretically,  based  on  a  conservation 
principle,  as  shown  in  the  work  of  Miralles-Wilhelm  et  al  (1989). 
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The  second  part  of  the  work  carries  the  study  of  the  mixing  phenomena  just 
described  into  situations  in  which  the  fluid  dynamics  are  more  complex.  Specifi 
cally,  the  intention  is  to  examine  the  process  of  molecular  mixing  when  it  occurs 
during  the  formation  of  coherent  structures  (vortex  rollup)  and  during  the  pairing 
and  merging  of  vortical  structures  as  it  takes  place,  for  example,  in  a  mixing  layer. 

The  inviscid  fluid  dynamics  of  vortex  pairing  and  merging  is  investigated  in  two 
different  situations  employing  the  method  of  vortex  dynamics  (Rangel  and  Si  rig- 
nano,  1989 ) .  The  first  situation  considers  the  pairing  and  merging  of  isolated  vor¬ 
tical  structures  and  of  vortical  structures  in  a  shearing  flow.  Thus,  the  mechanism 
of  pairing  is  analyzed  in  an  idealized  fashion.  The  second  situation  considers  vor¬ 
tex  rollup  and  pairing  of  two  or  grouping  of  more  of  these  structures  when  there  is 
an  initial  surface  of  discontinuity  or  vortex  sheet  separating  two  streams  traveling 
at  different  velocities.  The  temporal  instability  of  this  vortex  sheet  is  an  approxi 
mation  to  the  temporal  and  spacial  instability  occurring  in  a  turbulent  mixing  layer. 

The  method  of  solution  employed  in  this  study  is  based  on  the  dynamics  of  vor- 
ticity  motion.  Generally,  these  methods  take  advantage  of  the  fact  that  in  many 
flows,  vorticity  is  concentrated  in  regions  which  are  small  compared  to  the  total 
flow  domain.  The  evolution  of  the  vorticity  field  can  then  be  obtained  by  solving  a 
system  of  dynamic  equations  for  a  large  number  of  representative  vortex  cylinders. 
The  use  of  finite-core  vortex  cylinders  was  introduced  by  Chorin  and  Bernard  (1972) 
as  a  means  of  reducing  the  numerical  instability  inherent  to  this  problem. 

In  an  idealization,  the  mixing  layer  may  be  conceptualized  as  an  infinite  vortex 
sheet  subjected  to  a  certain  disturbance.  The  disturbance  grows  to  produce  the  finite 
vortical  structures  while  the  mechanism  of  viscous  diffusion  spreads  the  vorticity  on 
the  smaller  length  scale.  The  interaction  of  these  vortical  structures,  induced  by  the 
existence  of  subharmonics  of  the  initial  disturbance,  results  in  the  pairing  and  merg- 
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ing  of  t lie  structures.  In  the  inviscid  limit,  the  vortex  sheet  rolls  up  so  that  singular 
concentrations  of  vorticity  with  one  wavelength  spacing  are  formed.  If  subharmonics 
exist,  these  vortical  structures  interact  and  pairing  can  occur  (Ho  and  Huang.  1982). 

An  isolated  vortex  line  in  a  two-dimensional  inviscid  flow  produces  a  circular 
flow  pattern  where  the  circulation  along  any  closed  path  containing  the  vortex  line 
is  a  constant.  When  two  of  these  vortex  lines  are  located  a  certain  distance  apart, 
the  resulting  motion  is  such  that  the  two  vortex  lines  rotate  with  constant  velocity 
about  the  point  located  at  the  center  of  the  segment  joining  them  (Batchelor.  1970). 
The  separation  distance  remains  fixed  so  that  no  merging  occurs  regardless  of  the 
separation.  This  is  an  indication  that  such  a  concentration  of  vorticity  does  not  lead 
to  merging.  If  the  vorticity  is  not  initially  concentrated  along  a  line  but  confined 
inside  a  cylinder  of  finite  cross  section,  merging  is  possible  as  will  be  shown  below. 
The  initial  distribution  of  vorticity  over  a  circular  cross  section  can  be  conceptual 
ized  as  made  up  of  cylindrical  tubes  (rings  in  the  cross  sectional  area),  so  that  the 
total  vorticity  is  the  sum  of  the  contribution  of  all  rings.  The  interaction  of  two  such 
vortex  tubes,  one  from  each  vortex  cylinder  is  show  in  Fig.C.6.  Rotation  is  in  the 
clockwise  direction.  The  characteristic  length  used  in  the  nondimensionalization  is 
the  radius  of  the  ring,  while  the  characteristic  time  is  a  circulation  time  formed  with 
the  characteristic  length  and  the  circulation.  The  number  of  elements  employed  in 
the  calculation  is  entered  at  the  upper  right  corner  of  each  frame  as  n  and  all  the 
calculations  are  started  with  n=20.  For  the  case  illustrated  in  Fig.C.6,  the  ring 
separation  from  center  to  center  is  two  ring  diameters.  As  seen  in  this  figure,  this 
separation  is  small  enough  so  that  mutual  attraction  occurs  leading  to  merger  of  the 
two  rings.  For  sufficiently  large  separation  (more  than  four  diameters  from  center 
to  center),  each  ring  perceives  the  other  ring  essentially  as  a  line  vortex  and  neither 
attraction  nor  merging  occurs.  In  a  developing  mixing  or  shear  layer,  there  is  a 
spacial  evolution  of  the  vortical  structures  in  the  direction  of  the  flow.  At  a  given 
time,  the  downstream  structures  are  older  and  therefore  larger  than  the  upstream 
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ones. 


The  behavior  of  two  structures  of  different  size  is  presented  in  Fig.C.7  for  two 
structures  of  the  same  total  circulation  but  of  a  size  ratio  of  2/3.  This  is  an  idealiza¬ 
tion  of  the  effect  of  viscous  diffusion,  which  would  tend  to  spread  the  vorticity.  thus 
reducing  its  strength  but  conserving  the  total  circulation.  As  seen  in  Fig.C.7.  the 
smaller  structure  is  eventually  entrained  by  the  larger  one.  This  behavior  should  be 
of  significance  in  the  mixing  of  a  passive  scalar  or  two  reactive  species,  as  a  structure 
containing  one  of  the  reactants  may  be  trapped  inside  another  structure  containing 
the  other  reactant  thus  affecting  the  rate  of  reaction. 

In  a  shear  layer,  the  vortical  structures  may  be  initially  so  far  apart  that  they 
perceive  each  other  as  point  vortices.  Eventually  they  would  grow  as  a  result  of 
viscous  diffusion  but  this  mechanism  alone  is  not  sufficient  to  effect  a  fast  merging 
of  two  or  more  structures.  It  is  the  shearing  action  that  is  responsible  for  bring¬ 
ing  the  structures  closer  together.  To  visualize  this  effect  we  consider  two  vortical 
tubes  initially  separated  by  a  distance  of  four  tube  diameters.  Such  an  arrangement 
would  not  produce  pairing  and  merging  because,  as  explained  above,  the  structures 
are  sufficiently  far  apart.  In  a  shear  layer,  however,  these  two  structures  are  not 
alone.  Assuming  that  only  the  most  unstable  disturbance  is  present,  without  any 
subharmonics,  there  would  be  one  such  structure  every  wavelength  along  the  layer. 
In  this  case  the  wavelength  is  equal  to  the  separation  distance.  For  the  purpose 
of  this  visualization,  all  the  other  vortical  structures  are  represented  as  fixed  line 
vortices.  The  two  vortex  tubes  are  initially  displaced  a  small  amount  which  is  posi¬ 
tive  for  the  structure  on  the  left  and  negative  for  the  structure  on  the  right.  Figure 
C.8  shows  the  evolution  of  the  two  central  vortex  tubes  under  the  influence  of  the 
infinite  row.  Initially,  the  two  structures  move  slowly  towards  each  other.  After 
a  dimensionless  time  of  400,  the  close  range  interaction  leading  to  pairing  begins, 
resulting  in  the  merged  structure  shown  in  the  last  frame. 
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The  rollup,  pairing,  and  merging  of  vortical  structures  formed  during  the  evo¬ 
lution  of  an  inviscid  vortex  sheet  subjected  to  a  sinusoidal  disturbance  and  one  or 
more  of  its  subharmonics  has  also  been  investigated.  In  this  case,  the  surface  of 
discontinuity  is  discretized  into  a  large  number  of  vortex  lines.  A  sinusoidal  distur¬ 
bance  and  any  desired  subharmonics  are  then  applied  to  the  interface.  Acton  (197G) 
studied  the  evolution  of  the  inviscid  vorticity  field  using  vortex  dynamics.  Since  no 
rediscretization  or  insertion  of  new  vortices  was  performed,  the  structure  of  the  eddy 
was  lost  once  it  was  formed.  With  the  method  employed  here,  the  initial  resolution 
can  be  maintained  to  the  extent  of  the  calculation.  In  the  following  figures,  lengths 
have  been  made  dimensionless  with  the  disturbance  wavelength,  while  time  has  been 
made  dimensionless  with  a  residence  time  defined  as  the  ratio  of  the  wavelength  to 
the  magnitude  of  the  velocity  discontinuity  at  the  undisturbed  interface.  Figure 
C.9  illustrates  a  situation  in  which  a  sinusoidal  disturbance  and  its  first  subhar¬ 
monic,  both  of  dimensionless  amplitude  0.025.  are  simultaneously  imposed  on  the 
interface  at  t=0.  The  number  of  vortex  elements  n  is  indicated  at  each  time  frame. 
Figure  C.9  shows  how  the  presence  of  the  first  subharmonic  forces  the  pairing  of 
every  two  structures  in  the  manner  described  in  the  experiments  of  Ho  and  Huang. 
No  complete  merging  occurs  because  there  is  no  diffusion  mechanism  present  to 
spread  the  vorticity  away  from  the  interface.  Figure  C.10  shows  the  results  for  a 
case  in  which  the  second  subharmonic  is  present  in  addition  to  the  fundamental.  In 
this  case,  grouping  of  every  three  structures  occurs  and  most  importantly,  the  trans¬ 
verse  thickening  of  the  layer  occurs  more  rapidly  since  there  are  now  three  structures 
rolling  up  together  (Note  that  the  x  and  y  scales  in  Fig.  C.10  are  different  from  those 
of  Fig.  C.9).  The  evolution  of  the  vorticity  distribution  for  the  case  of  Fig.  C.10  is 
shown  in  Fig.  C.ll,  where  the  vortex  strength  (circulation  per  unit  area)  is  plotted 
against  a  coordinate  along  the  interface.  The  rollup  of  vorticity  into  three  struc¬ 
tures  is  apparent  in  this  figure  as  well  as  the  absence  of  viscous  diffusion  and  merger. 
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Figure  C.12  shows  the  amplitude  of  the  disturbance  as  a  function  of  time  for  the 
case  of  a  single  disturbance,  a  disturbance  plus  the  first  subharmonic,  and  a  distur¬ 
bance  plus  the  second  subharmonic.  The  exponential  growth  rate  predicted  by  the 
linear  theory  is  valid  until  a  dimensionless  time  of  0.8  approximately.  From  then 
on,  the  single-disturbance  grows  slowly  to  a  maximum  and  then  actually  decreases 
due  to  the  strong  shear  generated  by  the  resulting  vortex  row  and  the  absence  of 
viscous  diffusion.  The  cases  containing  the  first  and  second  subharmonics  continue 
to  grow  nearly  linearly  until  a  dimensionless  time  of  3.  From  this  point  on,  the  case 
containing  the  first  subharmonic  reaches  a  maximum  amplitude  and  then  decreases 
also  due  to  the  absence  of  merger  in  the  presence  of  strong  shear.  The  amplitude  of 
the  case  containing  the  second  subharmonic  continues  to  grow  linearly  although  at 
a  faster  rate  induced  by  the  three  structures  that  begin  to  roll  up  together. 


D.  VAPQRia.VC-PRQr.LET  CALCULATIONS 

The  detailed  analysis  of  a  vaporizing  droplet  in  the  convective  field  with  vari¬ 
able  thermo-physical  properties  has  been  examined  in  Chiang.  Raju  and  Sirignano 
( 1989).  We  present  here  some  of  the  results  of  our  recent  study.  The  base  case  study 
is  selected  as  a  cold  n-octane  fuel  droplet  suddenly  injected  into  the  hot  gas  stream. 
The  values  of  physical  parameters  in  the  base  case  are  given  in  Table  D-l.  Several 
parameter  studies  by  changing  initial  droplet  temperature,  ambient  temperature, 
initial  Reynolds  number,  and  droplet  heating  model  are  also  made  to  investigate 
the  effect  of  the  above  quantities. 

Typical  contour  plots  of  mass  fraction,  temperature,  vorticity  as  well  as  liquid- 
phase  streamlines  and  the  gas-phase  velocities  are  presented  in  Figure  D.  1.  The 
convective  effect  is  apparent  by  the  fore- aft  asymmetry  in  each  plot.  We  may  ob¬ 
serve  the  gas-phase  boundary  layers  surrounding  the  droplet  up  to  136°  in  angular 
direction  and  a  near  wake  region  behind  the  droplet.  The  vorticity  contours  are  con- 
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centrated  at  the  front  portion  of  the  droplet  and  then  are  convected  downstream. 
The  temperature  and  mass  fraction  contours  do  show  the  presence  of  the  boundary 
layer  and  wake.  However,  they  are  not  similar  to  each  other  because  of  non-unit  ary 
Lewis  number.  The  liquid-phase  streamlines  show  a  single  large  vortex,  roughly 
resembling  a  Hill’s  spherical  vortex,  in  the  interior  of  each  droplet.  The  similarity 
between  the  isotherms  and  streamlines  in  the  liquid  phase  shows  the  influence  of 
internal  circulation. 

The  surface  shear  stress  distribution  shown  in  Figure  D.  2. a  indicates  the  shear 
stress  decreases  with  time.  This  can  be  realized  by  the  following  facts:  (1)  With 
the  onset  of  surface  motion,  the  velocity  gradient  at  the  droplet  surface  decreases. 
(2)  The  surface  blowing  effect,  which  increases  the  thickness  of  the  boundary  layer 
and  reduces  the  velocity  gradient,  is  growing  as  vaporization  becomes  stronger.  (31 
As  droplet  heating  continues,  the  mass  fraction  at  the  surface  also  keeps  increasing, 
thus  yielding  lower  values  of  viscosity  of  the  mixture. 

The  tangential  velocity  distribution  shown  in  Figure  D.  2.b  seems  to  vary  with 
sinO  which  is  in  qualitative  agreement  with  Hill's  vortex  solution.  The  surface  tan¬ 
gential  velocity,  originally  zero,  is  brought  up  to  a  maximum  by  the  shear  stresses 
and  then  decreases  as  shear  stresses  diminish. 


The  distribution  of  vorticity  along  the  droplet  surface  is  very  similar  to  that  of 
shear  stress  as  shown  in  Figure  D.  2.c  The  surface  vorticity  is  obtained  from 
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The  dominating  term  is  $  which  is  also  the  main  term  contributing  to  the 
shear  stress.  It  is  clear  that  lj  ~  j  ~  Re ^2.  Hence  it  is  expected  that  the  effects 
of  internal  circulation,  boundary  layer  blowing  and  reduction  of  gas-phase  Reynolds 
number  will  cause  the  surface  vorticity  to  diminish.  If  the  flow  separation  point  is 
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determined  by  the  point  where  vorticity  changes  sign,  it  is  observed  that  the  sepa¬ 
ration  point  is  approximately  located  around  144  degrees. 

Figure  D.  2.c  shows  the  pressure  distribution.  The  recirculating  wake  dissipates 
part  of  the  kinetic  energy.  As  a  result  ,  the  pressure  cannot  recover  to  the  stagnation 
value.  We  now  switch  our  attention  to  the  overall  characteristics  of  momentum, 
heat  and  mass  transfer  to  the  droplet.  Figure  D.  3  shows  the  variation  of  drag 
coefficients  of  droplets  initially  at  three  different  ambient  temperatures  as  a  function 
of  hydrodynamic  diffusion  time.  Also  shown  in  this  figure  is  the  drag  correlation 
of  Renksizbulut  and  Haywood  (19%)  which  is  based  on  the  numerical  results  of  an 
isolated  moving  drop,  vaporizing  in  its  own  fuel-vapor  medium;  their  analysis  also 
considers  the  effect  of  variable  thermophysical  properties.  The  R-H  correlation  is 
given  by. 

CD(\  +  ^w./)°'2  =  ~— •(  1  +  0.2/fc^;63);  10  <  Rtm  <  300  (5) 

^  f  m 

where 

CH.;  =  CW(J~-r->(  ,_<?.) 

The  time  required  for  the  flowfield  to  relax  from  the  initially  impulsive  motion, 
due  to  the  sudden  injection  of  the  drop  into  a  uniform  flowfield.  can  be  estimated 
as  Adefax  ~  77^  ~  0.4  of  the  diffusion  time.  During  the  initial  relaxation  period 
the  drag  coefficient  falls  rapidly.  Subsequently  Cd  increases  as  a  result  of  a  reduc¬ 
tion  in  the  Reynolds  number.  However  for  the  higher  ambinent  temperature  case 
where  the  large  heat  flux  makes  droplet  vaporization  grow'  rapidly  during  the  early 
portion  of  life  time,  the  drag  coefficient  is  therefore  significantly  decreased  due  to 
the  boundary  layer  blowing.  For  all  cases,  Cd  tends  to  increase  with  time  during  the 
final  portion  of  the  calculation  which  implies  that  the  reduction  in  Reynolds  number 
takes  control  in  determining  drag  coefficients.  The  R-H  drag  correlation  seems  to  be 
not  applicable  in  our  cases,  except  at  the  early  times  of  the  lower  transfer  number 
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case  (ambient  temperature  =  800° A')  where  the  vaporization  is  weak.  The  deviation 
becomes  large  for  the  high  transfer  number  case. 

Three  components  of  drag  coefficients  are  shown  in  Figure  D.  3.  It  is  clear  that 
the  major  differences  in  total  drag  coefficients  for  three  different  ambient  temper¬ 
ature  cases  are  resulting  from  differences  in  the  friction  drag  coefficient.  A  higher 
transfer  number  can  reduce  the  friction  drag  by  a  large  amount.  The  thrust  drag  is 
always  negligible  and  at  most  accounts  for  less  than  2 %  reduction  of  the  total  drag. 
The  pressure  drag  coefficients  increase  steadily  which  may  be  due  to  the  increase  of 
the  wake's  strength  as  the  vaporization  increases  with  time. 

Tigure  D.  4. a  shows  the  average  A'usselt  numbers  of  these  three  cases  and  their 
corresponding  R-H  N’usselt  number  correlations  which  can  be  expressed  as 

Au/(  1  +  Bhj)°  51  =  2  +  OMRc'JiPr'/3 

The  Nusselt  number  falls  during  the  initial  relaxation  period.  Further  reduction 
in  the  Nusselt  number  is  believed  to  be  attributed  mainly  to  a  Reynolds  number 
reduction  and  an  increase  of  surface  blowing.  Again  our  numerical  results  do  not 
agree  with  the  R-H  correlation  as  in  the  drag  calculation.  However  the  discrep¬ 
ancy  becomes  smaller  during  the  final  calculation  when  the  surface  temperature 
approaches  the  wet-bulb  temperature.  The  discrepancy  in  both  C<j  and  A u  number 
between  our  calculation  and  R-H’s  calculation  may  be  attributed  to  the  different 
definition  of  the  problem  considered.  In  their  calculation.  R  and  H  considered  a 
droplet  vaporizing  in  its  own  fuel-vapor  environment.  As  a  result,  the  mass  diffu¬ 
sion  was  neglected.  In  addition,  they  keep  constant  droplet  surface  temperature  at 
the  boiling  point  which  may  reduce  the  heat  transfer  from  the  gas-phase  and  may 
cause  quite  different  thermophysical  properties  at  the  interface  when  compared  to 
our  calculation. 
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The  Sherwood  number,  which  represents  the  transient  dynamics  of  the  mass 
transfer,  is  shown  in  Figure  D.  4.b.  The  general  trend  of  variation  is  very  similar 
to  that  of  Nusselt  number.  During  the  early  relaxation  period,  the  large  mass  frac¬ 
tion  gradient  caused  by  the  suddent  injection  of  a  droplet  to  the  hot  gas-stream  is 
smeared  out  gradually.  This  relaxation  contributes  mainly  to  the  early  reduction  of 
Sherwood  number.  A  more  general  decrease  of  Sherwood  number  is  then  observed, 
which  may  be  due  to  the  decay  of  Reynolds  number. 

In  order  to  study  the  effect  of  initial  droplet  temperature,  the  initial  droplet 
temperature  has  been  raised  from  300/i  to  400/\.  The  results  are  shown  in  Figure 
D.  5. a  It  is  very  clear  that  for  the  higher  initial  droplet  temperature  case,  the  surface 
blowing  effect  is  more  significant  than  that  of  the  lower  droplet  temperature  case 
during  a  large  portion  of  the  droplet  life  time.  Hence  the  former  experiences  lower 
drag  than  the  latter.  Figure  D.  5.b  shows  the  portion  of  total  heat  flux  that  goes 
into  heating  of  the  drop  interior.  The  droplet  with  higher  initial  temperature  would 
spend  most  of  the  available  heat  energy  in  evaporation  process.  On  the  contrary,  the 
droplet  with  lower  initial  temperature  has  to  spend  most  of  the  available  energy  to 
heat  up  itself  first,  with  the  remaining  energy  going  for  vaporization.  As  a  result  the 
heating  process  is  slow  and  will  persist  for  most  of  the  lifetime.  The  vaporization  is 
also  relatively  weaker  compared  to  the  high  initial  temperature  droplet.  Also  shown 
in  this  figure  is  the  solid  sphere  correlation  which  can  be  grossly  inaccurate  due  to 
the  high  mass  transfer  at  the  surface. 

The  effect  of  variable  thermophysical  properties  is  always  interesting  for  the 
range  of  temperatures  considered  in  the  present  study.  A  calculation  has  been  con 
ducted  by  assuming  constant  properties  that  are  based  upon  the  values  in  the  far 
stream.  Since  the  global  behavior  is  evaluated  at  the  interface  where  the  mixture 
constituent  is  very  different  than  the  pure  gas  in  the  far  stream.  The  results  in  Fig¬ 
ure  D.  6  show  that  the  constant-property  calculations  result  in  the  overprediction  of 


the  drag  coefficient  by  about  20%.  The  discrepancy  is  attributed  to  the  change  in 
thermo-physical  properties  at  the  interface  and  the  change  in  the  flow  field  caused 
by  the  property  gradients. 

The  effect  of  initial  Reynolds  number  was  also  studied  during  the  current  work. 
The  calculations  for  three  different  initial  Reynolds  numbers  have  been  made.  The 
results  are  presented  in  Figure  D.  7.  It  is  observed  that  lower  initial  Reynolds  num¬ 
ber  results  in  a  higher  drag  coefficient  at  early  times  due  to  the  smaller  convective 
momentum  transport  which  usually  causes  the  early  separation  to  occur.  It  is  also 
noteworthy  that  the  three  drag  coefficients  at  the  final  calculation  tend  to  approach 
a  certain  form  which  is  governed  by  the  reduction  in  Reynolds  number  only. 

The  droplet  model  is  very  important  to  combustor  design.  Figure  D.  8  show 
the  comparison  results  from  an  infinite  conductivity  liquid-phase  model  where  the 
liquid-phase  temperature  is  assumed  to  be  uniform  in  space  but  varying  with  time. 
The  results  indicate  that  the  infinite  conductivity  case  predicts  higher  drag  in  the 
early  part  of  the  droplet  lifetime  and  lower  drag  later  on.  For  the  infinite  conductiv¬ 
ity  model,  the  liquid  phase  has  to  distribute  all  the  available  heat  flux  to  the  droplet 
interior.  Thus  the  surface  temperature  rises  slowly  resulting  in  very  weak  vaporiza¬ 
tion  rate  as  well  as  in  a  high  drag  coefficient  during  the  early  lifetime.  Because  of  the 
larger  driving  temperature  potential,  the  droplet  eventually  receives  more  heat  flux 
than  that  calculated  from  our  model.  As  a  result,  the  surface  temperature  builds  up 
and  the  vaporization  grows  quickly  while  the  drag  coefficient  is  reduced  significantly. 

With  the  results  of  above  cases,  it  is  desirable  to  obtain  the  correlations  of 
drag  coefficient,  Nusselt  numbers  as  the  function  of  some  important  parameters. 
Since  these  transport  quantities  are  determined  by  the  combined  influence  of  surface 
blowing  and  vaporization,  the  internal  circulation,  and  the  unsteady  effects  related 
to  the  droplet  deceleration.  A  complete  correlation  which  covers  a  wide  range  of 
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all  parameters  is  not  easy  to  achieve.  The  parameters  and  functional  form  selected 
in  this  work  are  the  same  as  in  the  R-H's  correlation.  The  following  modified 
correlations  show  a  good  agreement  (within  3%  discrepancy)  in  our  results  as  shown 
in  Figure  D.  9 

CD(  1  +  Bnjum  )0  336  =  -^-(1  +  0.225 Rc°m6' ) 

Nufllm(\  +  Bujiim )°  51  =  2  +  0.6  RcK2Pt){1 

for  .6  <  Buj.im  <  6/or)c/30  <  Rtm  <  250 

The  major  difference  between  the  our  modified  correlation  and  R-H  correlation 
is  the  exponent  of  ( 1  +  £//,/,(„,)•  Our  correlation  has  a  higher  exponent  which  em¬ 
phasizes  more  reduction  in  drag  coefficient  by  the  effect  of  transfer  number. 

The  droplet  interaction  effects  are  determined  by  the  combined  influence  of  sur¬ 
face  blowing  and  vaporization,  the  internal  circulation,  and  the  unsteady  effects 
related  to  the  droplet  deceleration.  Since  the  problem  contains  a  multidimensional 
variable  space,  an  exhaustive  numerical  study  of  all  possible  combinations  of  param¬ 
eters  would  be  extremely  elaborate.  Therefore,  we  have  concentrated  our  attention 
in  detailed  study  of  the  interaction  effects  rising  from  the  variation  of  Reynolds 
number  Reg_ o,  initial  diameter  ratio  a 2,o/aj,0’  an<^  initial  spacing  to  lead  droplet 
diameter  ratio  do-  Most  of  the  results  are  presented  in  Raju  and  Sirignano  (1987, 
1988).  The  main  observations  from  this  research  can  thus  be  summarized  as  follows. 

Droplet  collision  appears  to  be  likely  for  the  initially  equal-sized  droplets,  as 
long  as  the  initial  spacing  is  sufficiently  small  for  the  wake  effects  to  be  impor¬ 
tant.  However,  the  decrease  in  the  downstream  drop  size  can  prevent  collision  since 
droplet  deceleration  is  inversely  proportional  to  its  radus.  The  results  indicate  that 
there  exists  a  bifurcation  point  depending  upon  the  critial  ratio  of  a^o/^i.o  below 
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which  droplet  collision  becomes  unlikely;  the  smaller  the  ratio,  the  faster  the  rate 
at  which  drops  separate  front  each  other;  the  larger  the  ratio,  the  faster  the  rate 
at  which  drops  approach  each  other.  Closer  to  the  bifurcation  point,  the  droplets 
spend  more  time  under  the  influence  of  each  other  before  either  droplet  cllision  or 
separation  eventually  occurs.  Larger  initial  Reynolds  numbers  as  well  as  smaller 
initial  droplet  spacing  tend  to  reduce  the  critical  downstream  drop  size  required  for 
droplet  separation.  The  time  required  for  droplets  to  collide  is  relatively  insensitive 
to  the  Reynolds  number. 


The  lead  droplet  experiences  higher  drag  force  than  the  downstream  droplet 
which  is  shielded  by  the  lead  droplet.  The  lead  droplet  behaves  as  an  isolated 
vaporizing  drop  as  long  as  the  two  drops  are  separated  sufficiently  (2  <  d).  After 
studying  some  cases  the  drag  coefficients  of  the  two  droplets  can  be  correlated  by 


a2 


Cd, 2  =  “  l-22( 

V  oc ,2  fll  ^  oc  ,2 


)d0 


-0.25 


(6) 


The  numerical  results  satisfy  the  above  equation  within  8  percent  accuracy  for  the 
range  of  initial  spacings  between  2-15  and  initial  Reynolds  numbers  between  50-200. 


Figure  D.  10  shows  the  overall  Nusselt  numbers  Nu  and  average  surface  temper¬ 
ature  T,  for  the  case  of  two  approaching  droplets.  The  difference  in  both  Nu  and  T , 
between  two  droplets  indicate  that  wake  effects  tend  to  reduce  significantly  the  heat 
transport  to  the  downstream  drop.  For  the  case  of  separating  drops,  the  Nusselt 
number  for  the  downstream  drop  during  the  final  period  of  calculation  approaches 
the  theoretical  value  for  an  isolated  droplet  as  shown  in  Figure  D.  11. 

While  a  general  understanding  is  obtained  through  a  parametric  study  of  the 
interaction  effects  between  two  vaporizing  droplets,  the  results  are  by  no  means  con¬ 
clusive,  e.g.,  particle  dynamics  may  be  quite  different  if  a  very  small  particle  gets 
trapped  in  the  standing  eddy  of  the  lead  droplet.  Further  studies  are  needed  for  a 
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complete  understanding  of  the  drop  behavior  near  the  bifurcation  point .  The  present 
analysis  does  not  allow  for  droplet  coalescence,  or  atomization.  Future  studies  of 
the  dynamics  during  collision  would  be  worthwhile.  The  interactive  droplets  must 
also  be  studied  with  variable  properties  and  general  (three-dimension)  orientation. 
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Table  D-l.  Values  of  Physical  Parameters  Used  in  Base  Case  Computation 
for  Variable  Property  Single  Droplet 


Parameter 

Value 

Initial  Reynolds  number,  gas  phase  Re  =  toP'g  i00  /  p'oo 

100.0 

Relative  velocity  of  drop,  m/s 

25.0 

Free  stream  temperature,  °I< 

SOO.O 

Combustor  pressure,  atm 

10.0 

Prandtl  number,  gas  phase 

0.799 

Prandtl  number,  liquid  phase 

9.32 

Schmidt  number,  gas  phase 

2.23 

Molecular  weight,  oxidizer,  Kg/Kmol 

29.0 

Molecular  weight,  fuel,  n-octane,  Kg/Kmol 

114.2 

Droplet  initial  temperature,  °K 

300.0 

Viscosity  ratio,  p'J  p'goo 

9.30 

Density  ratio,  p',/p'goo 

362.47 

Specific  heat  at  constant  pressure  ratio, C’PJ  0/Cp'g  ^ 

1.69 

Latent  heat  /  Specific  heat  of  liquid,  °K 

136.8 
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FIGURE  A-l.  Plot  of  Distance  to  Breakup  versus  Shear  Air  Velocity  for 
Various  Liquid  Velocities. 
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FIGURE  A-2  Plot  of  Distance  to  Breakup  versus  Shear  Air  Velocity  for  One 
and  Two-Sided  Shear  Air  Flow. 


Figure  B.l.  Schematic  of  the  problem  of  liquid  sheet  atomization. 
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Figure  B.2.  Sheet  evolution  for  a  thickness-to-wavelength  ratio  of  2. 
Density  ratio=l.  Surface  Tension=0.  (shifted-symmetric  mode) 
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Figure  B.3.  Sheet  evolution  for  a  thickness-to- wavelength  ratio  of  1/2. 
Density  ratio=l.  Surface  Tension=0.  (shifted-symmetric  mode) 
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Figure  B.4.  Sheet  evolution  for  a  thickness-to- wavelength  ratio  of  1/4. 
Density  ratio=l.  Optimum  dimensionless  wavenumber,  (shifted-symmetric  mode) 
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Figure  B.5.  Sheet  evolution  for  a  thickness-to-wavelength  ratio  of  1/4. 
Density  ratio=l.  No  surface  tension,  (symmetric  mode) 


Figure  B.6.  Sheet  evolution  for  a  thickness-to-wavelength  ratio  of  1/8. 
Density  ratio*!.  No  surface  tension,  (symmetric  mode) 
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Fig.  C.l.  Schematic  of  the  problem. 
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Figure  C.6.  The  interaction  of  two  vortex  tubes  whose  centers 
are  initially  two  diameters  apart 
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Figure  C.7.  The  interaction  of  two  vortex  tubes  with  a  size  ratio  of  2/3. 
The  centers  are  initially  four  small  diameters  apart 


Figure  C.8.  The  interaction  of  two  vortex  tubes  whose  centers  are  initially  four 
diameters  apart  in  the  presence  of  an  infinite  row  of  point  vortices. 
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Figure  C.ll.  Evolution  of  the  vorticity  distribution  along  the  interface 
during  the  grouping  of  three  vortical  structures. 


MASS  FRACTION 
TIME  •  1.000 

CAS  PHASE  REYNOLDS  NUMBER  • 

AMIIIBT  TIimHATUM  -  Ml  hfV  U  » 
■mim.  »»0PUI  TtMPt«ATV»C  .  Ml 


GAS  AND  LIQUID  PHASE  ISOTHERMS 

TIME  •  >  000 

CAS  PHASE  REYNOLDS  NUMBER  •  SUIT 

•MtttaT  TIKHWTVM  •  m  hfM  Um 
BH1M1  Mcrui  tlurtUTmu  -  «h  dpi  uw» 


Figure  D.l.  Contour  plots:  (»)  pass  fraction,  (b)  temperature  of  liquid  and 
gas  phases,  (c)  stream  function  of  liquid  phase,  (d)  Isotherms  of 
liquid  phase,  (e)  vortlclty,  (f)  Instantaneous  velocity  of  gas 
phase. 
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Figure  D.2. 


(a)  Surface  shear  stress  distribution,  (b)  surface  tangential  velocity 
distribution  at  different  times,  1  -  time  -  0.98,  Re  «  91.30;  2  -  time 


-  3.00,  Re  -  93.83;  3  -  time  -  5.63,  Re  -  59.02;  4  -  time  -  8.30,  Re 
49.21;  5  -  time  -  11.25,  Re  -  38.61;  6  -  time  -  14.38,  Re  -  30.96;  9 
time  ■  19.13,  Re  -  25.63. 


Figure  D.2.  (c)  Surface  vorticlty  distribution,  (d)  surface  nondimensional  pressure 

distribution  at  different  times,  1  -  time  •  0.98,  Re  ■  91.30;  2  -  time 
-  3.00,  Re  -  93.83;  3  -  time  -  5.63,  Re  -  59.02;  4  -  time  -  8.30,  Re  • 
49.21;  5  -  time  -  11.25,  Re  -  38.61;  6  -  time  -  14.38,  Re  -  30.96;  9  - 
time  -  19.13,  Re  -  25.63. 
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Figure  D.3.  (a)  Time  variation  of  drag  coefficients  for  different  ambient 

temperatures,  (b)  time  variation  of  three  drag  components  for 
different  ambient  temperatures. 
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Figure  D.5. 


(a)  Time  variation  of  drag  coefficients  for  different  initial 
droplet  temperatures,  (b)  time  variation  of  the  rate  of  droplet 
heating  for  different  Initial  droplet  temperatures. 
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Figure  D.7.  Time  variation  of  drag  coefficients  for  different  initial 
Reynolds  numbers. 
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Figure  D.9.  (a)  Comparison  of  drag  coefficients  between  numerical  results  and 

correlations,  (b)  comparison  of  Nusselt  numbers  between  numerical 
results  and  correlations. 
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Figure  D.10.  Time  variation  of  average  Nusselt  number  and  droplet  surface 
temperature  for  approaching  droplets. 
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Figure  D.ll.  Time  variation  of  average  Nusselt  number  and  droplet  surface 
temperature  for  separating  droplets. 


AVERAGE  DROPLET  SURrflCE  TEMPERATURE 


